- Nom : Prénom : REAL Philippe
- Nom : Prénom : GUYONVARCH Alexis
Présentation des données et objectifs de l’étude https://www.kaggle.com/uciml/pima-indians-diabetes-database
Le jeu de données provient de l’Institut national du diabete, des maladies digestives et rénales. Il rend possible la prédiction de la pathologie, en l’espèce le diabète, pour le patient à partir d’analyses inclues dans le jeu de données. Dans l’extraction, les patients sont des femmes et issues de la commmunauté des indiens Pima.
Objectif : prédire si l’individu a ou non le diabète. Préalablement, il s’agi de compléter les données manquante par imputation.
## [1] "Pourcentage de données manquantes par variable"
## npreg glu bp skin Insuline bmi ped age
## 0.0 0.7 4.6 29.6 48.7 1.4 0.0 0.0
## type
## 0.0
Suppression de la 1ère colonne comprenant les identifants. La colonne insuline est, dans l’immédiat, conservée en dépit de la part importante de données manquantes : 48.7 %.
## [1] "Combinaisons 2 à 2 des données manquantes"
## glu bp skin Insuline bmi
## glu 5 0 0 4 0
## bp 0 35 33 35 7
## skin 0 33 227 227 9
## Insuline 4 35 227 374 10
## bmi 0 7 9 10 11
## [1] "Création des catégories de données manquantes"
##
## No-Missing Insuline.only bmi.only
## 392 140 1
## glu.only bmi.Insuline bp.Insuline
## 1 1 2
## bp.skin.Insuline bp.skin.Insuline.bmi glu.Insuline
## 26 7 4
## skin.Insuline skin.bmi.Insuline
## 192 2
Le pourcentage de données manquantes est élevé, les informations sont exhaustives pour seulement 51% des individus, ce qui jusitifie le recours à l’imputation multiple. Le graphique des combinaisons confirme, ce qui avait déjà été mis en lumière ci-devant, à savoir la prédominance de plusieurs combinaisons : * Insuline + Skin * Insuline + Skin + BP
A elle seule, la variable Insuline, quand elle est manquante, regroupe 8 patterns. La variable skin concerne 4 patterns.
Au final, il ressort que le mécanisme des données manquantes, qui concernent 5 variables du dataframe “Pima”, est non-monotone, ce qui justifiera ultérieurement le recours à l’imputation multiple (joint modeling, fully conditionnal specification, ACP).
La valeur d’ “influx” de la variable Insuline est plus élevée que la valeur d’“influx” de la variable Skin en dépit d’une proportion plus importante de données manquantes. Cela suggère une connection plus forte aux variables observées. Une valeur d’“outflux” très faible nous indique que la variable “Insuline”, ainsi que la variable “Skin”, quoique dans une moindre mesure, seront potentiellement moins utiles à l’imputation des autres variables.
##
## Variables sorted by number of missings:
## Variable Count
## Insuline 0.486979167
## skin 0.295572917
## bp 0.045572917
## bmi 0.014322917
## glu 0.006510417
## npreg 0.000000000
## ped 0.000000000
## age 0.000000000
## type 0.000000000
## missingCateg 0.000000000
## [1] "Nombre de patterns si Insuline manquant :8"
## [1] "Nombre de patterns si skin manquant :4"
On reprend la méthode d’analyse PCA des données manquantes vu en cours (chap 1) et de http://factominer.free.fr/missMDA/appendix_These_Audigier.pdf On ne remarque pas de groupes de variables bien séparés. Peut-être une liaison skin-insuline
graphiques N°1 - Mécanismes des variables “Insuline” et “skin” (§2.1)
Distributions marginales des variables Insuline, Skin et bp avec les autres variables incomplètes
Distributions marginales des mêmes variables avec les variables complètes
Dans un premier temps on traite les valeurs manquantes par imputation simple avec le package MICE : * par le biais de la méthode PMM (predictive mean matching), * puis au moyen d’une régression linéaire - non bayésienne - stochastique.
A priori, ces imputations ne seront pas conservées en raison de la sensibilité à la spécification du modèle (pour la méthode paramétrique) et du biais souvent généré par la méthode PMM (semi paramétrique).
#PMM
imp.si.pmm <- mice(data, m=1, seed = 111119, print = F)
#régression stochastique avec bootstrap
imp.si.norm <- mice(data, method = "norm.nob", m = 1, maxit = 1, print = F)
https://datascienceplus.com/imputing-missing-data-with-r-mice-package/
Nous recourons, au cours d’une première étape, à l’imputation multiple au moyen des méthodes Joint Modeling puis Fully Conditionnal Specification.
#JM
imp.mi.jm <- mice(data, m=5, method='norm',seed=111119, print=F)
#FCS
imp.mi.fcs <- mice(data, m=5, seed=111119, print = F)
Les graphiques montrent que les imputations par régression linéaire stochastique comportent des valeurs aberrantes, en l’espèce, des valeurs négatives pour l’épaisseur de peau et le taux d’insuline. Les distributions des valeurs observées et imputées sont proches, ce qui n’infirme ni ne confirme le mécanisme MAR. Pour les graphiques cf. [annexe N°1 - graphiques N°2 - Analyses des distributions pour les Imputations simples].
A nouveau, certaines imputations ne sont pas plausibles : avec le modèle joint les graphiques montrent des valeurs aberrantes pour l’épaisseur de peau et le taux d’insuline. S’agissant de la méthode FCS, les distributions des variables glu et bmi présentent des profils très divergents pour chacune des imputations, en raison du nombre peu élevé de données manquantes. Pour les graphiques cf. annexe N°1 - graphiques N°3 - Analyses des distributions pour les Imputations multiples.
Préalablement à l’étape de réduction des dimensions, nous procédons à un essai de transformation des variables visant à renforcer la linéarité des liens entre elles. Finalement, après plusieurs tentatives (log, logistic et racine carrée), seule une transformation logarithmique de la variable insuline est conservée.
Avec la méthode de validation croisée “kfold”, le nombre de dimensions optimale est de 6 ou 5 selon le critère “MSEP”.
## [1] "Le nombre de dimensions retenu pour le jeu de données initial est de :6"
## [1] "Le nombre de dimensions retenu pour le jeu de données transformé est de :5"
Visuellement, l’imputation de la variable insuline semble de meilleure qualité, une fois celle-ci ayant fait l’objet d’une transformation en logarithmique. Le nombre d’itération pour la période de burn-in Lstart et le paramètre L n’ont pas été déterminé par faute de temps en suivant la méthode proposé dans (http://factominer.free.fr/missMDA/appendix_These_Audigier.pdf) Or la détermination de ces paramètres joue un rôle important dans la qualité de l’algorithme Bayésien. Ceci explique peut-être le résultat moyen obtenu avec les paramètres par défaut.
Visuellement, le modèle semble plus en adéquation avec l’imputation multiple par bootstrap.
Cette dernière imputation semble in fine la plus appropriée. Plus de 90% des réimputations se situent dans l’intervalle de confiance.
L’analyse des distributions des résultats est mise en peuvre avec les méthodes du package “MICE”. L’analyse requiert une conversion préalable des objets avant manipulation des fonctions et méthodes du package. In fine, ces dernières sorties semblent confirmer le diagnostic précédent. Le meilleur modèle d’imputation semble être celui implémenté section 3.2.3.
L’analyse de sensibilité ne concerne que le modèle d’imputation multiple Avec pour hyptohèse un mécanisme MAR, les imputations sont robustes à une modification des valeurs de “skin” (idem pour “Insuline”). Seulement possible avec MICE (post traitement à la volée).
delta <- c(0,10,100)
imp.all <- vector("list", length(delta))
post <- mice(data[,1:8], maxit = 0)$post
for (i in 1:length(delta)){
d <- delta[i]
cmd <- paste("imp[[j]][,i] <- imp[[j]][,i] +", d)
post["skin"] <- cmd
imp <- mice(data[,1:8], post = post, maxit = 5, m=5, seed=111119, print = F)
imp.all[[i]] <- imp
}
bwplot(imp.all[[1]])
bwplot(imp.all[[2]])
bwplot(imp.all[[3]])
A partir des données complétées, la variable réponse, “Individus diabétiques”, est modélisée avec un modèle GLM logit avec pour prédicteurs : -npreg -glu -bp -skin -Insuline -bmi -ped -age
#Avec imputation MICE FCS
diabete.insuline.glm1 <- with(data=imp.mi.fcs ,exp=glm(type ~ npreg + glu + bp + skin + Insuline + bmi + ped + age,family=binomial(link="logit")))
#Avec imputations missMDA
diabete.insuline.glm2 <- with(data=conv.imp.mipca.boot.new ,exp=glm(type ~ npreg + glu + bp + skin + log.Insuline + bmi + ped + age,family=binomial(link="logit")))
diabete.insuline.glm <- with(data=conv.imp.mipca.EM.boot.new ,exp=glm(type ~ npreg + glu + bp + skin + log.Insuline + bmi + ped + age,family=binomial(link="logit")))
#Avec imputations missMDA mais sans la variable Insuline
diabete.noinsuline.glm <- with(data=conv.imp.mipca.EM.boot.new ,exp=glm(type ~ npreg + glu + bp + skin + bmi + ped + age,family=binomial(link="logit")))
summary(pool(diabete.insuline.glm2))
summary(pool(diabete.insuline.glm))
summary(pool(diabete.noinsuline.glm))
On remarque que les variables skin (pour mice), age et bp ont une p-value assez importante.
On accepte l’utilité de la variable “Insuline” et on conserve le modèle complet. p-value importante => on rejette l’influence de la variable insuline.
#Rapport de vraisemblance
pool.compare(diabete.insuline.glm, diabete.noinsuline.glm, method = "likelihood")$pvalue
## [1] 2.396047e-05
#Anova
anova(diabete.insuline.glm,diabete.noinsuline.glm)
## test statistic df1 df2 df.com p.value riv
## 1 ~~ 2 17.46443 1 374.3012 759 3.647477e-05 0.5507066
Au regard des différentes comparaisons et du test final, on conserve, tant pour la méthode d’IM FCS (package MICE) que celle utilisant l’algorithme EM et le bootstrap (package missMDA), des modèles excluant les variables : bp et age
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable age : 0.0763371217407416"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable skin : 0.017665142862979"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable insuline : 2.39604727777509e-05"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable bp : 0.0277648269328263"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable bmi : 0.0114624223656384"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable glu : 1.738609256563e-13"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable npreg : 3.1724938396982e-05"
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable ped : 0.00323902343173565"
Toutes les variables semblent significatives, c’est le modèle complet qui a le meilleur AIC. A partir de ce AIC moyen on pourrait essayer de supprimer les variables : bmi, age et bp.
## [1] "AIC du modele complet : 687.872740691428"
## [1] "AIC du modele sans la variable age : 689.177761785113"
## [1] "AIC du modele sans la variable skin : 695.753913386109"
## [1] "AIC du modele sans la variable insuline : 715.424610330563"
## [1] "AIC du modele sans la variable bp : 690.883997273214"
## [1] "AIC du modele sans la variable bmi : 694.39336092397"
## [1] "AIC du modele sans la variable npreg : 703.017985200401"
## [1] "AIC du modele sans la variable glu : 725.753952742626"
## [1] "AIC du modele sans la variable ped : 695.574620108385"
On regarde maintenant et l’on compare certains modèles emboîtés en supprimant successivement les variables : age, bmi et bp. Le meilleur modèle qui minimise le critère BIC est le modèle diabete.insuline.glm sans les variables : bmi, age et bp Soit en ne conservant que les 5 variables: npreg, ped, glu, skin et log.Insuline C’est aussi confirmé par les différents MLE réalisés sur les modèles emboîtés. Pour argumenter ce choix on peut se référer à l’annexe 2 - choix de modèle
## [1] "Test du rapport de vraisemblance modele complet vs modele sans la variable age : 0.0763371217407394"
## [1] "Test du rapport de vraisemblance modele sans la variable age vs modele sans la variable age et bp: 0.0860802995887788"
## [1] "Test du rapport de vraisemblance modele sans la variable age et bp vs modele sans la variable age, bp et bmi: 0.0424498515374343"
## [1] "AIC du modele complet : 687.872740691428"
## [1] "AIC du modele sans la variable age : 689.177761785113"
## [1] "AIC du modele sans la variable age et bp : 690.502897697801"
## [1] "AIC du modele sans la variable age, bp et bmi: 694.564436934921"
## [1] "BIC du modele complet : 729.666848289757"
## [1] "BIC du modele sans la variable age : 726.328079650294"
## [1] "BIC du modele sans la variable age et bp : 723.009425829835"
## [1] "BIC du modele sans la variable age, bp et bmi: 722.427175333807"
En reprenant (http://factominer.free.fr/missMDA/appendix_These_Audigier.pdf) (page 20) à propos de la sortie de la fonction pool. La colonne (fmi) (fraction of missing information) peut être interprété comme la part de variablilité due aux valeurs manquantes. De grande valeurs (>5) indique que le résultat est sensible à la méthode MI utilisée. Ici ce n’est pas le cas.
Class: mipo m = 100
estimate ubar b t dfcom df riv lambda fmi
(Intercept) -10.88502375 7.231029e-01 3.239346e-01 1.050277e+00 762 345.8641 0.45245827 0.31151206 0.31545909
npreg 0.15107510 8.032095e-04 5.978041e-05 8.635877e-04 762 683.0321 0.07517120 0.06991556 0.07262705
glu 0.02602838 1.685869e-05 2.779750e-06 1.966624e-05 762 574.4618 0.16653414 0.14275977 0.14572876
skin 0.04820190 7.002272e-05 3.368601e-05 1.040456e-04 762 329.4700 0.48588328 0.32699963 0.33104812
log.Insuline 0.90936111 2.642584e-02 1.944355e-02 4.606383e-02 762 242.1640 0.74313588 0.42632126 0.43100122
ped 1.00577056 8.573742e-02 5.498982e-03 9.129139e-02 762 695.2185 0.06477885 0.06083784 0.06352801
Au final, les estimations des coefficients sont du même ordre de grandeur en comparaison de la méthode des cas concrets (= listwise deletion). Mais de manière général sont plus faibles. Et en particulier biensûr pour les variables à fortes données manquantes: insuline et skin.
##
## Call:
## glm(formula = type ~ npreg + glu + skin + I(log(Insuline)) +
## ped, family = binomial(link = "logit"), data = data)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.0404 -0.6731 -0.3701 0.6458 2.5513
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -8.835302 1.183463 -7.466 8.29e-14 ***
## npreg 0.133916 0.039773 3.367 0.00076 ***
## glu 0.036197 0.005641 6.417 1.39e-10 ***
## skin 0.040632 0.013161 3.087 0.00202 **
## I(log(Insuline)) 0.250019 0.250644 0.998 0.31852
## ped 1.091993 0.400249 2.728 0.00637 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 498.90 on 392 degrees of freedom
## Residual deviance: 354.25 on 387 degrees of freedom
## (375 observations deleted due to missingness)
## AIC: 366.25
##
## Number of Fisher Scoring iterations: 5
Dans notre cas, nous avons réussi à démontrer la pertinence de l’imputation avec pour objectif l’implémentation d’un modèle de classification. Les méthodes sont nombreuses et nous n’avons pas pu toutes les expérimentés. Notament les méthodes bayésiennes de mipCA en ce qui concerne la détermination des paramètres. Ainsi que le package AMELIA couplé avec Zelig qui n’ont pas été utilisé. L’analyse comparative des données complétés par imputation multiple avec les packages mice et missMDA donnent des structures similaires, lorsque l’on fait une analyse PCA (voir Annexe 2).
Une question est apparue concernant l’utilisation pratique des données imputés. De nombresue méthodes R nécessite en input un datafRame ou équivalent. Ici les méthode utilsant du boostraping, on a pour résultat une famille de structure comparable à un dataFrame. Que faut-il faire pour se ramener à un dataFrame: utilser la moyenne? Ou par exemple dans le cas la méthode mipca utiliser la variable $res.InmputPCA c’est ce qui est fait en Annexe 2 - pour utilser les méthodes de choix de modèles et en annexe §3 pour une analyse PCA des données complétées.
Les données Pima sont très utilisées (cf. Kaggle). Pour autant il ne semble pas qu’un modèle d’inputation ait été utilisé. Bien souvent les données manquantes sont supprimées (analyse des cas concrets ou listwise deletion) ou bien une imputation simple par moyenne est mise en oeuvre.
Distributions marginales des variables Insuline, Skin et bp avec les autres variables incomplètes
Distributions marginales des mêmes variables avec les variables complètes
Les graphiques montrent que les imputations par régression linéaire stochastique comportent des valeurs aberrantes, en l’espèce, des valeurs négatives pour l’épaisseur de peau et le taux d’insuline. Les distributions des valeurs observées et imputées sont proches, ce qui n’infirme ni ne confirme le mécanisme MAR.
Imputations simples - PMM
Imputations simples - Normales
A nouveau, certaines imputations ne sont pas plausibles : avec le modèle joint les graphiques montrent des valeurs aberrantes pour l’épaisseur de peau et le taux d’insuline. S’agissant de la méthode FCS, les distributions des variables glu et bmi présentent des profils très divergents pour chacune des imputations, en raison du nombre peu élevé de données manquantes.
Imputations multiples - JM
Imputations multiples - FCS
Dans cette annexe on applique les techniques regsubset et step pour le choix automatique des variables. Les données utilsées sont les données obtenu à partir de la variable $res.inputePCA On obtient à nouveau ce résultat en utilsant comme données d’input les moyennes des diffrérents jeux de données issues du boostrap résultant de mipca.
La méthode step du package leaps nous fait choisir un modèle à 4 variables explicatives: npreg + glu + skin + ped
On reprend le modèle saturé obtenu au §2: m_sature
m_sature = glm(formula = type ~ . , family = binomial(link="logit"), data = res.imputePCA[,1:maxDim])
summary(m_sature)
##
## Call:
## glm(formula = type ~ ., family = binomial(link = "logit"), data = res.imputePCA[,
## 1:maxDim])
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.3533 -0.6207 -0.3008 0.6030 3.0587
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -13.582666 1.151711 -11.793 < 2e-16 ***
## npreg 0.151973 0.035413 4.291 1.78e-05 ***
## glu 0.014893 0.004842 3.076 0.00210 **
## bp -0.018386 0.008602 -2.137 0.03256 *
## skin 0.056456 0.014345 3.936 8.30e-05 ***
## log.Insuline 1.623134 0.227638 7.130 1.00e-12 ***
## bmi 0.028698 0.022177 1.294 0.19565
## ped 0.983273 0.304624 3.228 0.00125 **
## age 0.021588 0.010604 2.036 0.04175 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 993.48 on 767 degrees of freedom
## Residual deviance: 633.74 on 759 degrees of freedom
## AIC: 651.74
##
## Number of Fisher Scoring iterations: 5
Critère AIC
Critère BIC
#modele_step_BwdFwd_AIC<-glm(formula = type ~ (npreg + glu + bp + skin + Insuline + bmi + ped + age), family = binomial, data = dataFramePima)
summary(modele_step_BwdFwd_AIC)
##
## Call:
## glm(formula = type ~ npreg + glu + bp + skin + log.Insuline +
## ped + age, family = binomial(link = "logit"), data = res.imputePCA[,
## 1:maxDim])
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.3924 -0.6149 -0.3075 0.6151 3.1102
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -13.304451 1.124008 -11.837 < 2e-16 ***
## npreg 0.145590 0.034766 4.188 2.82e-05 ***
## glu 0.015047 0.004842 3.107 0.00189 **
## bp -0.015841 0.008386 -1.889 0.05888 .
## skin 0.069347 0.010370 6.687 2.27e-11 ***
## log.Insuline 1.646173 0.226420 7.270 3.58e-13 ***
## ped 1.018496 0.302159 3.371 0.00075 ***
## age 0.021184 0.010588 2.001 0.04542 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 993.48 on 767 degrees of freedom
## Residual deviance: 635.41 on 760 degrees of freedom
## AIC: 651.41
##
## Number of Fisher Scoring iterations: 5
AIC(modele_step_BwdFwd_AIC)
## [1] 651.4134
BIC(modele_step_BwdFwd_AIC)
## [1] 688.5637
#m_StepBwdFwd_BIC<-glm(formula = type ~ (npreg + glu + skin + Insuline + ped),family = binomial, data = dataFramePima)
summary(m_StepBwdFwd_BIC)
##
## Call:
## glm(formula = type ~ npreg + glu + skin + log.Insuline + ped,
## family = binomial(link = "logit"), data = res.imputePCA[,
## 1:maxDim])
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.3949 -0.6142 -0.3125 0.6117 3.0078
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -13.699976 1.040001 -13.173 < 2e-16 ***
## npreg 0.175578 0.029819 5.888 3.91e-09 ***
## glu 0.015125 0.004724 3.202 0.001364 **
## skin 0.062768 0.009899 6.341 2.28e-10 ***
## log.Insuline 1.654175 0.225200 7.345 2.05e-13 ***
## ped 1.027049 0.298013 3.446 0.000568 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 993.48 on 767 degrees of freedom
## Residual deviance: 641.26 on 762 degrees of freedom
## AIC: 653.26
##
## Number of Fisher Scoring iterations: 5
AIC(m_StepBwdFwd_BIC)
## [1] 653.2596
BIC(m_StepBwdFwd_BIC)
## [1] 681.1223
PCA sur les données complétées par MIPCA
PCA sur les données imputé par mice
PCA sur les données Pima de MASS (non imputées)
Le nuage des individus est bien centré
Dans le cas de mice toutes données sont relativement centrées.
Remarque: Les données complétés par MICE et MIPCA ont une PCA similaire. Par contre pour les données de MASS on a une transformation par rapport aux axes. Même forme mais symétrique / à l’axe principal.
A la vu des graphiques, la structure globale des corrélations est très voisine.